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SUMMARY 



Recently the fiber optic technology has been introduced to the petroleum industry and 
one of its main advantages is to get continuous measurements of temperatures in space 
and time inside the wells. There are different applications of temperature measurements; 
one of its main applications inside the producing wells is zonal rates allocation ftom 
different producing zones. Due to the weak relationship between the temperature and the 
zonal rates, so the main objective of this study is to develop a rigorous mathematical 
forward model of temperature as function of zonal rates and some other fluid and thermal 
properties starting from the basic balance equations. Then, use this forward model as a 
tool in inversion to invert the measured temperature to allocate the rates from the 
different producing zones. A simple forward model suitable for handy calculations has 
been obtained for two-phase (oil-water) two zones producing wells. The accuracy of the 
inversion for the zonal allocation from temperature measurements has been examined 
using synthetic examples for single and two-phase oil-water flow under different 
conditions. In addition the forward model has been tested on field data taken from a 
producing well where a fiber optic distributed temperature sensors (DTS) has been 
installed to test tile forward modeling and to validate the findings from the synthetic 
examples. The forward modeling shows a good agreement with the field data of 
maximum difference less than 1 F and the results of the inversion from the synthetic and 
field examples reveal the important of imposing the total rate as constraint during the 
inversion especially for low temperature contrast between the producing zones. 
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Chapter 1 
INTRODUCTION 

[1.1] Objectives 

Since the introduction of the permanent down hole fiber-optic sensors to the oil and gas 
industry in early 1990 and many works have been done to try to implement this 
technology and study its application. Since the last two years, Schlumberger has been 
deeply involved in the application of this technology 1 and they developed SENSA 
distributed temperature sensing (DTS) where the temperature measurements inside 
production and injection wells can be done continuously in space and time, temperature 
can be obtained every 7 seconds and collected every 3.3 ft 1 . These continuous 
measurements of temperature have many applications in the oil and gas industry, one of 
these applications is to know when and from where the production occurs in a producing 
wells and the most important is how much flow is taken place from multi-laterals 
producing wells. Accordingly, the objective of this work is to study the availability of 
using the DTS measurements in quantifying the zonal rates from multi-zones producing 
wells and test how accurate is these quantifications. 

[1.2] Working Steps 

To achieve the above objectives, the following steps have been taken: 

1- Developing a simple analytical forward modeling that describes the temperature 
profile in a single phase oil/water, single layer producing well that capture most of 
the physics of the problem 

2- Make a comparison with the other known models and with the numerical solution 
available in Eclipse 300 under the thermal option. 

3- Reformulate the forward model to be applied for a single liquid phase two layers 
production well and extend the modeling to two-phase oil and water production. 

4- Show the concept of the inverse modeling to invert the temperature measurements 
to estimate the rates from the different producing zones by the aid of some 
synthetic examples. The objective of these synthetic examples is to: 
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i- Study the availability to invert the temperature for the zonal rates 

under different conditions, 
u- Do uncertainty analysis to see the effect of the uncertainty of the 

model parameters on the accuracy of rate allocation from different 

zones. 

5- Test the application of the forward and the inverse modeling on a field case using 
DTS measurements in one of highly deviated production well 

So the organization of the report will be as follows, the second chapter will concern 
with the first two points, the third chapter will handle the third point, the fourth 
chapter will cover the fourth point, the fifth chapter will be mainly for the field case 
application, and finally, the sixth chapter will be devoted for conclusions and 
recommendation for future works, mainly the extension of the forward modeling to 
handle multi-phase production. 
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Chapter 2 



TEMPERATURE FORWARD MODELING IN SINGLE 
LAYER PRODUCTION WELLS & MODEL SELECTION 



The objective of this chapter is to show the mathematical formulation of the temperature 
modeling in vertical or deviated producing well from a single layer-single phase liquid 
production, propose different models by either including or neglecting some model 
parameters, and finally compare all the proposed models and the well known models 
from the literature with the numerical modeling obtained from Eclipse 300 under the 
thermal option in order to select the best model for a single layer production and then use 
this model for extending the modeling to more than one layer production. 

[2.1] Mathematical Formulation of the Temperature Modeling 

Fig. 2.1 shows the thermal nodal analysis used to develop a mathematical temperature 
model by determining the temperature at each node using mass, momentum and energy 
balance equations. 



Zyft 




Fig 2.1: Schematic diagram showing the nodal temperature analysis in tubing and in 
formation. 
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Table 2.1 shows the temperature nomenclature at each node presented in Fig. 2.1 . 



Table 2.1: Temperature Nomenclature at each node 



Node 


Nomenclature 


l 


Bottom noie tormation lemperaiure caicuiaieu trom me carm gcuuici uuu 
gradient, Tcibh 


2 


Bottom hole flowing fluid temperature, T ft hi 


2' 


Flowing fluid temperature at the top of the producing zone, Tfoh2 


3 


Formation temperature at the well/earth interface in the non-producing zone 


4 


Formation temperature calculating from the earth geothermai gradient in the non 
producing zone, Tei 


5 


Flowing temperature, Tf, and it can be at any depth Z from the producing zone j 



Depending upon the transport phenomena, whether mass, momentum, or energy, between 
each node, the balance equation is used. The following diagram, Fig. 2.2 shows the type 
of balance equation and the assumptions used in each equation. These assumptions are 
taken based on the previous works of many authors and also to try to keep the equation 
easy to implement and at the same time does not violate the physics of the problem. 
The material balance equation in general form is written as follows: 
Mass Balance Equation: 

Rate of increase of mass = rate of mass in - rate of mass out (2.0) 

Momentum Balance Equation: 

Rate of increase of momentum = rate of momentum in - rate of momentum out + external 

force on the fluid •••• (2M) 

Energy Balance Equation: 

Rate of change of (internal energy + K.E+P.E. due to convection)+(net rate of heat 
addition by conduction)-(net rate of work done by the system on the surrounding) = (Rate 

of accumulation of internal energy +KJE + P.E) ~(2.000) 

For the reader convenience, I highly recommend to return to see the general mathematical 
equation for the above three balance equations in Bird et al book 2 (Table 1 1.4-1, page 
341) 
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Analytical Modeling of Temperature Profile in a 
liquid phase (oil/water) pr duction wells 



d 

Producing Zone 
(Node 1-2) 



Energy 
Balance 



Assumptions: 

2- In production, temperature at 
the perforations (Tdbh) is the 
same as the temperature of the 
earth calculated from the 
geothermal gradient 

3- Conduction heat transfer is 
neglected 

4- Work done by the fluid 
against the viscous force is 
neglected 

5- Steady State Problem (No 
energy accumulation in the 
system) 

6- Change in P.E. is neglected 

7- Incompressible fluid and 
neglect the area change 
between the two nodes, so 
change in K.E is neglected 



Equation 2.8 



In front of the 
producing zone 
(Node 2 -D 



Non Producing Zone 
(Node 4 -3) 



I 



Energy + Mass 
Balance 



Energy Balance 



Assumptions: (Mixing 


between streams) 


1- 


Steady state (No 




mass or energy 




accumulation) 


2- 


Neglecting change 




inP.E.andICE. 


3- 


No loss or gain of 




heat during mixing 




(adiabatic mixing) 


4- 


Fluid is 




incompressible or 




compressibility is very 




small 


5- 


Work done by the 




fluid against the viscous 




force is neglected 


6- 


Mixing takes place 




at constant pressure 


7. 





Equations 2A 1& 2.12 



Fig 2.2: Analytical modeling of temperature profile in liquid 
production wells 



i. 



Well Path 
(Node? -51 

— T~ 



Energy + Mass + 
Momentum Balance 



Assumptions: 

1- Work done by the fluid against 
the viscous forces is neglected 

2- Thermal conductivity is 
constant 

3- Heat conducted from the 
producing zone is neglected 



i 



Equations 2.19 & 2.20 



1- Steady State Problem (No 
energy, mass and 
momentum accumulation in 
the system) 

2- Work done by the fluid 
against the viscous forces is 
neglected 

3- Constant heat flux form the 
tubing to the casing and 
from the casing to the 
surroundings at each control 
volume 

4- Thermal resistance of pipe 
and steel is neglected 
compared to that of the fluid 
in the tubing/casing annulus 

S* Incompressible fluid and no 
area change, so change in 
K.E. is neglected 



Equations 2.43, 2.38 
&2.39 
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[2.1.1] Producing zone (node 1 and 2) 

The general energy balance equation written in term of enthalpy instead of internal as 
presented in Bird et al 2 is as following: 

Rate of change of (enthalpy* KE.+PJZ due to convection)+(net rate of heat addition by 
conduction)-(net rate of work done by the system on the surrounding) = (Rate of 

accumulation of enthalpy + KB + PE) (2.1) 

According to the assumptions mentioned in Fig, 2.2 for node (1-2), the general energy 
balance equation, Eq. 2.1 will be: 

dH = 0.0 (2.2) 

From the basic thermodynamic principles, dH can be obtained from the following 
equation 3 , 

dH = C p dT-MjrC p dP ( 2 - 3 ) 

Where, fl^ is the Joule Thomson coefficient and can be obtained from the following 
equation 3,4 : 

'-[SR[#H 

For incompressible or slightly compressible fluid, p = cons tan t , V = cons tan t 

:0.0 (2.5) 



From Eqs. 2.5, 2.4, and 2.2, Eq. 2.3 will be 

dH = C p dT+VdP = 0.0 (2.6) 

:.C„dT = -— (2.7) 

P 

Solving Eq. 2.7 leads to Eq 2.8, which is written in the field unit as follows: 

144(P,-/V) 
P f C p J 

Eq. 2.8 indicates that only the effect of Joule Thompson coefficient is dominant in 
calculating the temperature at node 2 by knowing the temperature at node 1 due to the 
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flow through the perforations. The pressure drop can be calculated using Darcy's 
equation by knowing reservoir rock and fluid properties and assuming a steady state 
flow 5 : 



[2.1.21 In ^nt of the Producing zone (node 2 - 2*) 

In this part, to obtain an expression for the temperature, we need to divide the producing 
zone into equal intervals, each interval producing equal rate. The number of divisions 
depends upon the number of temperature measurements in the producing zone. By 
applying a macroscopic mass and energy balance 2 due to the mixing of two streams, we 
can get the temperature at the any interval inside the producing zone using the following 
derived equation: 

T 9m CpCM/Tfthfrj) +<li C p(o T d n 10) 

Where, i = 2, 3, , n (n is the number of temperature measurements inside the 

producing zone), taking into consideration that T^hi is calculated from Eq. 2.8 and 2.9 
before. Also, T* should be corrected due to the pressure drop across the perforation using 
the same Eqs. 2.8 and 2.9 but using Td instead of Teu>h 

As the fluids produced from each interval inside the producing zone have equal rate and 
equal specific heat capacity, Cp, so Eq. 2.10 can be written in the following form: 

t^ J-^™* 1 * (2 - n) 

Accordingly, temperature at node T will be: 

T Wfl =r jM(0) (2-12) 

It should be noted that Eq. 2. 1 1 is rate independent as it depends upon assuming that at 
each interval inside the producing zone, the producing rates are equal and the sum of 
those individual rates is the total producing rate from this producing zone, accordingly we 
should impose a condition such that at no production or physically at neglected 
production, Eqs. 2.10 and 2.1 1 does not hold and in this case the temperatures inside Ihe 
producing zone should be equal to the geothermal temperature. 
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For the synthetic cases, and to simplify the calculation, we can assume that only one 
temperature measurement inside the producing zone (one single entry) and in this case 
temperature at node T is the same as at node 2. 



[2.13] Non Producing zone (node 4 -3) 

As the fluid produced, heat is transferred by convection inside the well bore and some of 
this heat is lost by conduction to the non-producing formation. Thus, inside the non- 
producing zone, the transport phenomenon is only the heat energy due to heat loss from 
the well bore to the non-producing zone. So the only balance equation required is the 
energy balance equation. 

By applying the general energy balance equation given in Eq. 2.1 between node 3 and 4 
and according to the assumptions mentioned in Fig. 2.2 for node (4-3), the equation 
becomes: 

ar 1 \ ar ~ k s dt {2 U) 

It should be mentioned that Eq. 2.13 is in ID radial with the most common assumption 
that the earth density is unvaried with space and also with constant earth thermal 
conductivity. 

Eq.2.13 can be converted to the dimensionless form by using the following dimensionless 



Where Th is the temperature at node 3 

r D =— (2.15) 

<»=777r<=7-< ("9 

Pe C c r «6 r "* 

The radial partial differential equation of temperature distribution in earth in 
dimensionless form will be: 
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9 2 r eD | 1 37^ ^37^ £.17) 

Eq.2.17 can be solved using the following initial and boundary conditions: 

Initial Condition: limT^ =0, Temperature is constant (equal Td, which is the earth 

tp-*0 

temperature at any given depth and at infinite distance away from the well, temperature at 
node 4). 

Boundary Conditions: 

dT 

Outer Boundary Condition: lim — — = 0 (No change in temperature at infinite r D ) 



dT 

Inner Boundary Condition: — & 



= -1 (rate of flow of heat from the well bore to 
the surrounding earth across an element dz is constant). 

Eq. 2.17 is converted to the Laplace domain and is solved using Mathematica software. 
The solution of T eD in the Laplace domain (Teo(s)) is in the following form: 

Where, s is a dummy variable for the Laplace domain, and lo, Ii, Ko, Ki are the modified 
Bessel functions. Eq. 2.18 is the cylindrical source solution of Eq. 2.17 which is difficult 
to be inverted to the time domain analytically, so I used the Gaver functional, Wynn-Rho 
algorithm that is coded in Mathematica to get the inversion numerically by setting r e o to 
a very high value (e.g. 1000) and r D to 1 as we are interested to get the T c d at the 
well/earth interface. Fig. 2.3 shows the Tcd versus t D at the well/earth interface. 
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0.00001 

1/10- 10 l.'l(r 7 0.0001 0.1 100 100000. 

t D .dimensionless 

Fig. 2.3: Log-Log plot of Ted versus t D at the well/earth interface 

Hasan & Kabir 6 have done an approximation to the solution of Eq. 2.17 at the well/earth 
interface using the same initial and boundary condition and got the following equation: 
1^ =1.1281^1.1-0.3^] if t D £ 1.5 (2.19) 

r eD | ro= ,=[0.4063+0.51n(r o )]jl+^j ift D >1.5 (2-20) 

A comparison between the solution using the numerical laplace inversion and that 
obtained from Hasan & Kabir 6 solution is given in Fig. 2.4, the comparison shows that 
most of the points range from 10' 10 days to 10000 days lies on a 45° line, so for 
simplicity, Hasan and Kabir 6 solution, Eqs. 2.19 and 220 is used to get the temperature at 
node 3 by knowing the temperature at node 4 obtained from the geothermal gradient. 
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0 



2 3 4 5 6 7 
Ted solution using numerical laplace inversion 
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Fig. 2.4: Comparison between the numerical laplace inversion solution and Hasan & 
Kabir solution for Ted at the well/earth interface 

[2.1.4] Well Path (node T- 5) 

As the fluid proceeds from node T to 5, heat energy is transported by convection and 
also mass and momentum are transported due to the fluid flow. So, energy, mass, and 
momentum balance equations are applied between node 2* and 5. 
The general energy and mass balance equation after applying the assumptions mentioned 
in Fig. 2.2 will be: 

^Lj£_ v ^- g Sin* (2.21) 

dz dz dz 

For radial heat transfer from the well bore fluid to the well (cement) /earth interface, 



Where, T h is the temperature at node 3 and w t is the total mass flow rate, which is 
calculated as follows 3 : 



dQ_ ln r a U 
dz w, 



(T f -T h ) 



.(2.22) 
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' ~ 1.1309xl0 6 



246.6 



.(2.22a) 



The radial heat transfer from the well (cement)/earth interface to the surrounding can be 
obtained from Eq. 2.14 which is the definition of die dimensionless earth temperature, 



dz 



.(2.23) 



.(2.24) 



From Eq. 2.22 and 2.23, 

dQ Jln [ r s UK c 1 _ 
dz w t \K,+T*r t u\ Kf 

From the basic thermodynamic principles or from Eq. 2.3, the specific enthalpy can be 
taken from the following formula: 



dH 



dp 



dT, 



.(2.25) 



By substituting Eq.2.25 and Eq. 2.24 in Eq. 2.21, the final energy balance equation for 
the fluid in the producing well will be: 

Cp £L = _*£_ [ r t UK. l T dp. dv 

" dz w, [K e +7^ r„ U] K 1 " 5 ™ p dz dz 

Where, U is the overall heat transfer coefficient and can be calculated from Eq. 2.27 as 
shown in (Bird et al) 2 under the following conditions according to Sagar et al 3 : 

a- Thermal resistance of pipe and steel are negligible compared to the thermal 
resistance of the fluid in the tubing/casing annulus, 

b- Radiation and convection coefficients are negligible and can be ignored 



U = 



.(2.27) 



Eq. 2.26 is considered as a general equation for the temperature distribution inside the 
well bore between node 2' and 5 after combining both mass and energy balance equation 
and according to the assumptions mentioned in fig. 2,2. Eq. 226 has been given by 
different authors 3, 7 and can be applied for both single and multi-phase flow. 
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By applying the momentum and mass balance equation according to the assumptions 



mentioned in Fig. 2.2, the term (^) can be obtained as follows: 

dz 



dp dv 

— =-0V 

dz 



pv- — pgsin0-^H 
dz dz 



•(2.28) 



Substituting Eq. 2.28 in Eq. 2.26, 



(2.29) 

The analytical solution to Eq. 2.29 depends upon the joule Thomson coefficient, which is 
difficult to be obtained for a multi-phase flow. So, to make, the problem easy, we have to 
divide the solution to a single phase liquid (oil or water production), two phase liquid (oil 
+ water production), single-phase gas, and multi phase (oil + water + gas). This chapter is , 
mainly concern to give the analytical solution for single-phase liquid production and the 
other chapters will give the analytical solution for the two-phase liquid production (oil + 
water). 

For the single phase liquid, I can consider black oil production below the bubble 
point pressure especially because I am considering only on die analysis of DTS very 
close to the producing intervals in which the gas hold up is usually very small 
compared to that on the surface, so that the assumption of constant density and that 
tile pressure is below the bubble point pressure is applicable at a very small interval 
close to the producing zone for a black oil production. From Eq. 2.4 for 

V = cons tan /, the joule Thomson coefficient becomes: 
-1 



-v- 



dv 
dz 



Mjt = 



PC, 



,(2.30) 



By substituting Eq. 2.30 into Eq. 2.29, and write the equation in field units, 



dT. 



dz 



la 



w,C p 



r„UK e 



144 dp 



pC. dz 



....(2.31) 
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The pressure loss due to friction can be obtained from the equations mentioned in Beggs 
and Brill 8 or the technical description of Eclipse 300 9 as follows: 
2.956x10-" fpq 2 



dz 



Where the friction loss coefficient, £ can be obtained as follows: 



If if e <2000, /=— , ifU e >2000 -L=-3.61og 



R, 13.7DJ 

Where Re is the Reynolds number and can be obtained as follows: 
_ 0.1231pq 

By substitute Eqs. 2.32, 2.33, and 2.34 into Eq. 2.31, Eq. 2.31 will be: 



dT f In 



dz w,C p 



K * +T "n u ] 



Where, ^H4. 2.956x10" fq^ 

T* is calculated by knowing the earth temperature at the bottom hole (Ten*), which is a 
fixed quantity, and the earth temperature gradient using the following equation: 
T^T^-GrzSme (2.37) 



.(2.32) 



..(2.33) 



.(2.34) 



..(2.35) 



.(2.36) 



Eq. 2.35 can be converted to a dimensionless form using the following dimensionless 
parameters: 



_ V 



Z D=~ 



Z 

L 



.(2.38) 
(2.39) 



4d- 



-In L 


r 6 UK t 


w,c p 





1.12x 86400 J 



.(2.40) 
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(2.41) 



By substituting Eq. 2.37 and the dimensionless foims, Eqs. 2.38 to 2.41, Eq. 2.35 will be: 

^•*.<r*-*j** a '"? z ' h M > < 242 > 

<^D I Jbh l jbh 

The boundary condition used to solve the above ordinary differential equation is 
Tpizjy = 0) = 1 This means that the temperature at the bottom hole is equal to Tfot,. This 
boundary condition is suitable for dealing with a production from a single layer. Another 
boundary condition should be used if dealing with the solution of Eq. 2.42 in multi layers 
producing wells, as we will see in the next chapter. It should be noted here that the origin 
of the z scale is at the bottom hole as shown in Fig.2.1. The solution of Eq. 2.42 is done 
using Mathematica to get the profile of the dimensionless fluid temperature inside the 
well in front of the non producing zones as function of the dimensionless depth and the 
solution is as follows: 

-G r SinJ L-A D G T Sinfl z D L^A D T^ h -r^^-hexp^z D )>[g r Sing L+A D {T M -T m )+T^ D ] 
Tp - A D \ h 

(2.43) 

Eqs. 2.38 and 2.39 are used to convert the profile from the dimensionless domain to the 
real domain by knowing the fixed fluid temperature at the bottom hole of the well, T^h, 
and the depth of the well, L. 

[2.2] Model Selection 

The major objective of this part is to compare the model proposed from the first part of 
this chapter with the existed work in the literature and also with the numerical solution 
obtained from Eclipse 300 under the thermal option to test the accuracy of the 
mathematical model proposed. Also, as we can see that the high uncertainness in the 
proposed model is during modeling the temperature in front of the producing zone 
(conrecting the geothermal temperature due to the Joule Thomson effect because of the 
pressure drop across the perforation). As in modeling the temperature in this part, we 
need reservoir information about the flow regime, and some reservoir properties like 
permeability, drainage radius and the producing thickness, where there is a huge 
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uncertainty in obtaining those parameters. So, the second objective of this part is to 
neglect the correction of the geothermal temperature in front of the producing zone due to 
the pressure drop in other words, neglect the joule Thomson effect in front of the 
producing zone and compare the results with the numerical and see how much the error 
will be. 

The models used for the comparison with the numerical solution using Eclipse 300 under 
the thermal option are as follows: 

1- Ramey 10 solution for the liquid phase production, 

2- Sagar model 3 , and 

3- The proposed model given by Eq. 2.43, the proposed model is named as "model 3" 
during the comparison. 

As also the second objective is to study the effect of including the Joule Thomson 
coefficient in front of the producing zone in the temperature modeling, so additional 
models have been included in the comparison, which are mainly obtained from Eq. 2.43. 
These models are: 

4- The proposed model, which is "model 3" by neglecting the Joule Thomson coefficient 
in front of the producing zone, this model is named "Model 4", 

5- Model obtained from Eq. 2.43 by neglecting the friction loss ( fa 52 0) and neglecting 
the Joule Thomson coefficient in front of the producing zone, this model named as 
"model 2", 

6- Model obtained from Eq. 2.43 by neglecting the friction loss ( = 0) and including 
the Joule Thomson coefficient in front of the producing zone, this model is named 
"model 1" 

In addition to the above six models, two other models have been obtained from Ramey 
and Sagar model by applying the Joule Thomson coefficient on those models in front of 
the producing zones, these models are named as follows: 

7- "Modified Ramey model* 

8- "Modified Sagar model" 

Table 2.2 shows the well and the reservoir fluid and rock data used by the above eight 
different models and the numerical model obtained from Eclipse 300 under the thermal 
option to do the comparison and select the best model 
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Table2.2: Input Data for the numerical (Eclipse 300) and the forward models: 



Well Data 


Value 


Units 


q a , oil Rate 


1020 


STB/D 


Ow. water rate 


0 


STB/D 


q Q , gas rate 


306 


MSCF/D 


GLR I 


300 


SCF/STB ! 


Well Length. L 


6792 


ft 


r H (inside tubing radius) 


0.9075 


in 


rto (outside tubing radius) 


1.1875 


in 


r d (inside casing radius) 


2.506 


in 


to, (outside casing radius) 


2.75 


in 


Tvah (well bore radius) 


3.75 


in 


Heat transfer coefficient between the well bore and the formation, 
U 


156.55 


BTU/D-fP-F 


Well roughness 


0.001 


ft 


oil heat capacity, Cpo 


0.465 


BTU/lbm-F 


oil API 


58 




Inclination to the horizontal, 6 


90 










Rock thermal Properties 






Thermal rock conductivity (Ke) 


33.6 


BTU/D-ft-F 


Heat capacity of rock 


0.21156 


BTU/lbm-F 


Geothermal Gradient (G T ) 


0.0274 


Fffl 


Temperature of the earth at the bottom hole of the well (Tam,) 


237.2 


F 








Reservoir Properties 






Permeability (homogenous reservoir), K 


200 


md 


thickness of the producing layer, h 


20 


ft 


drainage radius, r s 


0.62 


ft I 


oil formation volume factor, B„ 


0.91 


BBL/STB 


Simulation Time, t 


168 


hr 



Figure 2.5 shows the comparison of all the eight models with respect to the simulation 
results, all the models shows good agreement except for Sagar and Modified Sagar 
models due to the correlation used by Sagar 3 which might be not working for this case as 
the comparison was with Eclipse results not with a real data as this correlation was 
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developed from a real data or it might be the data used in this example are beyond the 
data base from which this correlation was developed 

Figure 2.6 shows the point-by-point comparison with Eclipse results for all the models 
except Sagar and modified Sagsr as they both show high error compared to the other 
models. This figure shows the plot of the error between each model and simulation 
results versus depth, as it is seen from figure 2.6 that model 3 and 4 give less error 
compared to the other models less than 0.3 F difference from the simulation results for 
the distance close to the producing layer about 600 ft (at depth of 6192 ft). 
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Fig. 2.5: Comparison of the eight analytical forward models with the numerical model 
(Eclipse 300 [ thermal option) 
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Fig.2.6: Error Comparison of the all the models except Sagar and Modified Sagar 



Above the 6192 ft depth, models 3 and 4 gave higher error compared to the other models 
about 7 F difference from the simulation results at the surface and about a maximum of 
0.6 F difference from the model of least error. The main reason that model 3 and 4 gave 
good results near the producing interval compared to the simulation is that gas holdup is 
still very small so the assumption of single phase liquid is holding, also the assumption of 
constant friction loss is reasonable whenever the gas holdup is very small which is the 
case near the producing interval. While at far distance from the producing interval the gas 
holdup is increasing and the assumption of constant friction loss and single phase liquid 
production does not hold anymore that is why model 3 and 4 are not giving good results 
at higher distance from the producing interval. 

As we are interested only on the temperature analysis near the producing zones as there is 
a high chance to have a high liquid holdup and the assumption of liquid phase can be 
achieved so model 3 and 4 are the best candidates. Since as seen from Fig. 2.6 that the 
difference between model 3 and 4 are small and could be neglected, so model 4 is 
selected in order to reduce the uncertainty of imperfect knowledge about the reservoir 
properties like the type of the flow regime, permeability, drainage radius, reservoir 
thickness, etc 
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Chapter 3 

TEMPERATURE FORWARD MODELING IN TWO 
LAYERS PRODUCTION WELLS 

The objective of this chapter is to extend the temperature forward modeling derived in the 
previous chapter to be applied for two or multi-layers production wells for both single 
and two-phase liquid production. So, the first part of this chapter will concern on the way 
of extending model 4, which was derived from Eq. 2.43 by neglecting the Joule Thomson 
coefficient in front of the producing zone to be applied for the two layers production and 
the second part is to give the way of extending the forward model to two phase liquid two 
layers production wells. 

13.1J Temperature Forward Modeling for single-phase liquid, two layers production 
wells 

Figure 3.1 shows a thermal nodal analysis sketch for two layers, single-phase liquid 
production wells. 




Fig. 3.1: Nodal Analysis for two layers, single-phase liquid production wells 
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The only difference between the single layer and the two layers production is in the nodal 
analysis between nodes 5-5\ nodes 8-7, and nodes 5*-9, also there is a minor change 
between nodes 2'-5. The main differences between the single and the two layers 
production will be mentioned for each of these nodes. 



[3.1.11 Node 2*-5 

The ordinary differential equation given across these nodes for the single layer 
production, Eq. 2.42, is the same as that used for the two layers production however the 
boundary condition that will be used is more general such as: T© (zd = Zdbh ) = Tabh 
This general boundary condition allow us to handle the two or multi layers production 
cases as the temperature at node 5* should be corrected due to the mixing between the 
two streams and also due to the change of the rate from qi to qi+q2, so in this case as if 
we are dealing with the well as consisting of different sections each has the same 
equation but different boundary condition depending upon the temperature of the 
previous section. 

The solution of Eq. 2,42 using the above general boundary condition has been done using 
Mathematica and the solution is as follows: 

{A D T^ -G r Sin0(A D z D L+L)-^ 0 D )+ 
* y*** (~A D T eibh +A D T tt T^+Gr Sin0(A D z m L+L)+T M 

(3.1) 

Where, Tfdbh is the temperature of entry and z^h is the depth measured from the bottom of 
the well at the entry level Eqs. 2.38 and 2.39 are used to convert the dimensionless 
temperature profile obtained from Eq, 3.1 to the real domain. 



r - i 

1° AT 

A D A fl>h 



[3.1.2] Node 5 -y 

The modeling between node 5 and 5' is very similar to that between node 2 and T for the 
single layer production presenting in the previous chapter in that both mass and energy 
balance are applied, also the assumptions used between node 2 and 2' are the same as 
between node 5 and 5* except the last assumption where the heat capacity of the two 
streams are not the same and also the mixing rates are not equal. 
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Similarly, by dividing the producing zone into equal intervals, each interval producing 
equal rate. The number of divisions depends upon the number of temperature 
measurements in the producing zone and by applying a macroscopic mass and energy 
balance 2 due to the mixing of two streams from the upper producing zone and the total 
rate obtained from the lower zone, we can get the temperature at any interval inside the 
producing zone using the following derived equation: 

T m =1 ?U (3-2) 

n n 

Where, i = 1, 2, n (n is the number of divisions or the number of temperature 

measurements inside the upper producing zone) 

Tfl): is the temperature at each interval inside the producing zone, Tffl» is the well bore 
temperature at node 5. 

Cpi: is the specific heat capacity of the fluid in the upper producing zone, 
qi, q2 : is the total production from the lower zone and upper zone respectively 
Cp(i): is the specific heat capacity at each interval inside the producing zone and is 
calculating as a rate weighting average according to the following equation: 

c _L 5J " (3.3) 

*°~ [.♦**] 

Where, i = 1, 2 ,n and C p( o)= Cpi which is the specific heat capacity of the fluid 

produced from the lower zone, which is constant through the section between node T and 
5. 

At node 5', 

T — T (3.4) 

Also, <V-C^ < 3 - 5 > 

It should be noted that unlike Eq. 2.11, which is used to model the temperature from the 
lower producing interval (node 2-2'), Eq. 32 and 3.3 are rate dependent, however, it 
should be mentioned also that as the total rate from the two producing zones are null or 
very close to zero, Eq. 3.2 and 3.3 will blow up, so a condition should be imposed such 
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that at very small or zero rates from the two producing zones, the temperature should be 
. equivalent to the geothermal temperature so that to honor the physics of the problem. 
Also, it should be noted that Td in Eq. 3.2 could be corrected due to the pressure drop 
across the perforation in a similar way mat has been described in the previous chapter by 
using Eqs. 2.8 and 2.9 or the Joule Thomson effect can be neglected and Td will be 
directly the geothermal temperature as was found from the previous chapter that the 
correction due to the Joule Thomson effect in front of the producing zone could be 
neglected. 

In case of synthetic examples, and to simplify the calculation, we can consider that we 
have only one temperature measurement at the top of the upper producing interval so that 
i =1 and n = 1 in Eqs. 3.2 and 3.3 so that we can get T 5 - by knowing T s and T 6 

13.13] Node 8-7 

The same equation described in the previous chapter between node 3 and 4 can be used 
between node 8 and 7 but the flow rate is the total rate from the two producing zones. 

[3.U1 Nodes'- 9 

Eq. 3.1 can be used to describe the temperature profile between node 5' arid 9 by using 
the total rate (qi + qz) instead of qi. Also, Cp between node 5* and 9 is equal to Cps- as 
calculated from Eq. 3.5. 
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(3.21 Temperature Forward Modeling for two-phase liquid, two layers production 
wells 

Fig. 3.2 shows the nodal analysis sketch for two phase two layers production, where the 
two phases are oil-water (liquid), so that the forward model derived before will be 
applicable. 



Non Producing 




Producing zone 2 



Producing zone 



Fig.3.2: Thermal Nodal analysis for two phase, two layers production well 



The extension of the modeling to two-phase flow depends upon recalculating the 
parameters of the modeling for the two-phase flow. The equation for each parameter will 
differ depend upon the nodal location, thus, the equation for each parameters will be 
given between each node with a special reference to the equation used in the temperature 
modeling. 

It should be noted that in the non producing zone as there is no fluid flow, only heat 
energy flow, so the change from single phase to two phase flow will not affect the 
temperature modeling between nodes 3 and 4 and node 8 and 7. 
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Also, it should be mentioned that the correction of the temperature due to the pressure 
drop in front of the producing interval is neglected otherwise we need to have and 
additional inputs to the model which are the relative permeability and the capillary 
pressure which will increase the uncertainty of the modeling. 

13.2.11 Node 2-2* 

As seen from Eq. 2.11, that the temperature modeling between these two nodes depends 
only on the geothermal temperature, which does not depend upon the production phase, 
so the temperature modeling between node 2 and T is the same as for single-phase flow. 



13.2.2] Node 2'- 5 

Eq. 3.1 is used to get the temperature distribution between these two nodes. The 
1 parameters with its equation that are required to be obtained due to the two-phase flow 
are as follows: 

1- For A D calculation the parameters required are: w t and Cp, w t is calculated using Eq. > 
2.22a where q w is substituted by q*i and q 0 is substituted by q„i, while Cp between 
nodes 2* and 5 is calculated according to the following equation: 

c _ C j» fr>i+ C P» fl»' (3.6) 

2- For the <J>d parameters the parameters that are required to be obtained for the two- 
phase flow are: q, Cp, p, p. Those parameters can be obtained as follows: 

i~q = q l =q ol +q M ( 37 ) 

ii- Cp is calculated as given by Eq. 3.6 

m-u^iL= Mo ' q ' i+Mw ' 9 * (3-8) 

j*pmg l m K m **+ p -'** (39) 

[3.2.31 Node 5-5* 

Eq. 3.2 and 3.3 are used to calculate the temperature between these nodes by substituting 
qi by (q 0 i + q w i ) and q* by (q 0 2 + qw2), Cp 2 is calculated from the following equation 
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c C^'q^+C^-q^ 

" 2 ~ q o2 +q„2 

[3.2.41 N de5*-9 

Same as between node T- 5 mentioned before but with the following differences: 

i- For w b it is calculated using Eq. 2.22a, where q w = q W T = qwi + q*2 and 

qo=qoT=qoi + qoi 

ii- For Cp, it is calculated using the following equation 

Hi- For the <fc> calculation the parameters required to be calculated are: 

ff-flrt+f* * (312) 

A, • (g a i + q,i) ± Mw • (g*i + .(3.13) 

?,l+?»l+9o2 + «w2 

= P. -fai+g.a)+fl> .(3.14) 



The extension of the temperature modeling to multi-layers two-phase flow is trivial as 
it is only an extension of the equations mentioned in this chapter and in the previous 
one. 
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Chapter 4 

TEMPERATURE INVERSE MODELING WITH 
SYNTHETIC EXAMPLES 

As the major objective of this study is to use the temperature forwarding modeling 
discussed in the previous two chapters as the forward tool and use it in inverting the 
temperature measurements inside the producing wells to allocate the rates from the 
producing layers. So, the main objective of this chapter is to give a short introduction 
about how we can do the inversion of temperature for rate allocation, test the 
applicability to invert the temperature for zonal allocation for a single phase production 
under different scenarios in order to pick the critical scenarios through which the 
inversion might or might not work. Then use these critical scenarios to study the accuracy 
of zonal allocation by introducing errors in the forward model parameters, and finally, 
testing the accuracy of the inversion for two phase production on the critical scenarios 
proposed before. So the organization of this chapter will be, first give a brief introduction 
about the inversion, second testing the inversion on different synthetic examples for 
single-phase two layers production, then show the accuracy of the inversion for the errors 
in the forward model parameters and finally testing the inversion for two phase two 
layers production wells. 

14.1] Inverse Temperature Modeling 

Inversion in its broad meaning is to find the independent parameters in the forward model 
that minimize the error between the measured dependent parameter and the calculated 
dependent parameter from the forward model. So, it is an optimization problem where we 
want to minimize a certain objective function, which is the error between the measured 
and the calculated dependent parameters, by changing the independent parameters in a 
certain domain, or in other words changing the independent parameters according to 
specified constraints. 

In this study, the dependent parameter is the temperature and the independent parameters 
are mainly the zonal rates or could be any other input parameters of the forward 
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modeling. Accordingly, the .mathematical description of the optimization problem is as 
follows: 



Where, m: is a vector of the independent parameters, mainly the zonal rates and some 
other input parameters. 

There are many optimization algorithms in different books 11, 12 that can be used to find 
the zonal rate by minimizing the error between the temperature measured from DTS 
(Distributed temperature sensors) and the calculated temperature from the forward 
modeling discussed in the previous two chapters, however for simplicity and to make the 
analysis of DTS fast and easy to handle for the users, I used the optimization algorithm 
coded in Excel which is the "Generalized Reduced Gradient" 13 so that the user can do the 
analysis using Excel without any programming skill 
[4.2] Inversion Testing for Single-Phase Two Layers Producing Wells 
The major concern of this part is to see whether the rate can be allocated by knowing 
only the temperature measurements in producing wells under different scenarios, the 
selection of these scenarios depends upon an early work 14 that shows that the main 
factors which affect the rate allocation are 1- temperature contrast between the two zones, 
2- the production time, and 3- the amount of the flow rates. Accordingly eight different 
scenarios have been tested using synthetic examples for a single phase two layer 
producing well to see which is the most severe factor that has high effect on the rate 
allocation. Those scenarios are as follows: 

1- Small rate from the lower zone (200 STB/D) and high rate from the upper zone 
(800 STB/D), low temperature contrast between the two producing zones (about 3 
F, equivalent to a separation distance of 100 ft), and low production time (10 
days). 

2- Small rate from the lower zone (200 STB/D) and high rate from the upper zone 
(800 STB/D), high temperature contrast between the two producing zones (about 
13 F, equivalent to a separation distance of 500 ft) and low production time (10 
days). 




(4.1) 



28 



3- Small rate from the lower zone (200 STB/D) and high rate from the upper zone 
(800 STB/D), low temperature contrast between the two producing zones (3 F, 
equivalent to a separation distance of 100 ft) and high production time (100 days). 

4- Small rate from the lower zone (200 STB/D) and high rate from the upper zone 
(800 STB/D), high temperature contrast between the two producing zones (13 F, 
equivalent to a separation distance of S00 ft) and high production time (100 days). 

5- High rate from the lower zone (800 STB/D) and low rate from the upper zone 
(200 STB/D), low temperature contrast between the two producing zones (3 F, 
equivalent to a separation distance of 100 ft) and low production time (10 days). 

6- High rate from the lower zone (800 STB/D) and low rate from the upper zone 
(200 STB/D), high temperature contrast between the two producing zones (13 F, 
equivalent to a separation distance of 500 ft) and low production time (10 days). 

7- High rate from the lower zone (800 STB/D) and low rate from the upper zone 
(200 STB/D), low temperature contrast between the two producing zones (3 F, 
equivalent to a separation distance of 100 ft) and high production time (100 days). 

8- High rate from the lower zone (800 STB/D) and low rate from the upper zone 
(200 STB/D), high temperature contrast between the two producing zones (13 F, 
equivalent to a separation distance of 500 ft) and high production time (100 days). 

Through all the above 8 scenarios, I used the input data given in Table 4.1 in addition to 
the rates, time and the information about the separation distance between the two zones to 
get a true temperature profile inside the producing well using the forward modeled 
developed in the previous two chapters. Then, starting with initial guess of zero or very 
small rates from the two producing zones, which is equivalent to assuming the initial 
guess is the geothermal temperature and assuming a perfect knowledge in the input data 
by using the same data set given in Table 4.1 in addition to the production time and the 
separation distance between the producing zones given in each scenario to do the 
inversion on only the true generated temperatures very close to the two producing 
intervals. So, the in Eq. 4.1 is the true generated temperature and T^ is the 
calculated temperature from the forward model using the same input data that is used to 
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generate the true temperature but with different rates, which are very small rates for the 
initial guess and there is no constraint imposed during working in those scenarios. 

Table 4.1: Input data for the synthetic examples of the different eight scenarios for single 



phase two layers production well 





Value 


Unite* 

units 


Formation Thermal Properties 






Geothermal gradient, Gt 


0.02 1 4 


r/u 


rhermal diffusivity. a 


0.04 


u £/nr 


Thermal conductivity, Kq 


oo.o 


oTt 1/n.ft-P 


Bottom hole temperature, T^h 
1 


237.2 


c 
r 


Fluid properties 






Oil API 


oo 


API 


Oil specific heat capacity, Cp 


0.485 


BTU/Jbm-F 


Oil viscosity, po 


1.06744 


cp 


Well Data 






Total length of the well, L (Depth of the lower zone) 


6792 


ft 


Inside radius of the tubing, r a 


0.9075 


in 


Outside radius of tubing, r to 


1.1875 


in 


Inside radius of the casing, r d 


2.506 


in 


Outside radius of the casing, fa 


2.75 


in 


Well bore radius, r^ 


3.75 


in 


Annulus water thermal conductivity, 


9.192 


BTU/D-ft-F 


Cementing material thermal conductivity, kcem 


96.5 


BTU/D-ft-F 


Tubing roughness 


0.001 


ft 


[Inclination anale to the horizontal, 6 


90 


Deg. 



All the results from all the scenarios show the high ability to invert easily for allocating 
the rates from the two producing zones by knowing only the temperature profile very 
close to the producing intervals. Figure 4.1 shows the results from scenario 7, as it is 
found that this scenario took the largest number of iterations to converge, about 14 
iterations using the generalized reduced gradient (GRG) coded in Excel and using only 
the points near the producing zones. In Fig. 4.1, the initial temperature profile is the 
temperature calculated from the forward model using very small rates as initial guess, the 
true temperature (considered as measured temperature) is the temperature generated from 
the forward 



30 



Tf(F) 






Initial 


Final 


TRUE 


qo1, STB/D 


0.001 


800.1194 


800 


qo2,STB/D 


0.001 


200.0314 


200 



Fig. 4 J: Scenario 7 results for rate allocation 

model using the true rates and the other data mentioned in each scenario and the final 
temperature is the temperature at final iteration, also given in the Fig. 4.1 the initial, true 
and final rates from each zone. It should be noted that in the above synthetic examples for 
each scenario model 4 is used as the forward model, in which the Joule Thomson 
coefficient is neglected in front of the producing zones only. There are two things that 
have been done for the synthetic examples for those 8 scenarios and also for all the 
synthetic cases that will be presented in this chapter just to simplify the calculations. 
First, there is no temperature points in front of the lower producing zone, so that node 2 
and T are considered on node 2\ second, there is only one temperature point in front of 
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the upper producing zone and its location is at the top of the upper zone, which is at node 
5\ 

From the above testing of the inversion for the different 8 scenarios, it was found that the 
only scenarios that show a big difference in the number of iterations are the scenarios that 
contain the high and low temperature contrast compared to the minor difference due to 
the time and rate changes. So, in the next part I will study the accuracy of the rate 
allocation by introducing error in the model parameters for the scenarios that have high 
and low temperature contrast. The scenarios selected for this study are scenarios 7 and 8 
as they both show higher number of iterations to converge compared to the others. 

[4.3] Testing the Accuracy of ?tate Allocation for Single-Phase Two Layers 
Producing Wells 

The major concern of this part is to show how the uncertainties in the forward model 
parameters affect the rate allocation from the two producing zones. Nine cases have been 
studied to show how the effect of the errors in the parameters affecting the rate allocation 
from the two producing zones. The forward model used is Model 4 and the effect of each 
case has been studied for scenarios 7 and 8. Table 4.2 and Fig. 4.2 show the different 
cases that have been studied to test the accuracy of the rate allocation for the low and 
high temperature contrast scenarios. 

It should be noted that the input data used to generate the true temperature for the nine 
cases are the same as that given in Table 4.1 and scenarios 7 and 8. So, case one, where 
we assume a perfect knowledge of the input data, is exactly the same as the results 
obtained from scenarios 7 and 8, while in the other cases, error has been introduced to 
most of the model parameters to see its effect in the accuracy of the zonal allocation. 
Case 9 shows a combined error of some of the parameters. Also, it is important to 
mention that mat the objective function used in those nine cases are the same as in Eq. 
4.1 using only the points near the producing zones and without constraint However for 
cases 8* and 9* the total rate has been introduced as constraint in the optimization 
problem and also using the points very close to the producing zones. 
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Table 4.2: Cases description and errors in rate allocation due to errors in model 
parameters for low and high temperature contrast 





Conditions 1 


Error in ol. % 
Low temp. 1 High temp. 
Contrast Contrast 


Error In < 
Low temp. 


»2. % j 
High temp. 1 


Comments j 


1 1 1 


Perfect information 1 
about all the | 
parameters 1 


0.0027 


0.00103 


ft AfVTA 
U.UU /** 






2 


200% error in the I 
roughness of the 
tubing I 


4,51 


2.67 


4.14 


2.8 j 


High error is used for the 1 
roughness as it is the only 
parameters that has a high 
uncertainty J 


j 3 1 


10% error in oil heat 1 
capacity j 


4.68 | 


9.15 


4.66 


9.05 




1 A 

1 4 j 


10% error in oil 
specific gravity j 


0.19 


1.67 


0.06 


1.66 1 






10%errorinthe 


0.5 1 


0.5 


0.7 


0.7 1 




6 


10% error in the 
overall heat transfer 
coetrtcient 


2.94 


4.08 


2.93 


4.06 




1 * 


10%errorinthe ' 
geothermal gradient 


8.3 1 


9.2 


3.7 


3.8 




8 


Error in temperature 

measurements 
(normal distribution 
of mean zero and 
1 SJ> of 0.1), in 
addition using only 
two digits precision 
1 for the temperature 


74.5 


8.18 


74 9 


15.18 




8* 


Case 8, including 
total rate as 
constraint 


5 


0.7 


20J 


3 


By including the total rate as 
constraint in the objective 1 
1 function, the results of the 1 
I error in the rate allocation 1 
has been improved | 


9 


Mixing of all the 


1 69.3 

r ' 


| 1.11 


72.4 


18 




9* 


Case 9, including 
total rate as 
constraint 


2.45 


3.7 


9.8 


14.8 


By including the total rate as 
constraint in the objective 
1 function, the results of the 
J error in me rate allocation 
] has been improved 
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Error Analysis 
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Fig. 4.2: Error analysis of rate allocation for error in the model parameters 



From the above results it was found that: 

1 - In all the cases that have error of 10 % in its parameters, the rate allocation from 
the two producing zones show error less than 10 % except for cases 8 and 9. Case 
8, where errors in the temperature measurements has been generated which ranges 
from -0.3 to 0.3 F, while case 9 including all the error in the parameters in one 
single case 

2- Improvement in the rate allocation has been shown after including total rate as 
constraint in the objective function. Also, another observation is that case 9, 
which includes all the error in the parameters shows approximately good results 
compared to case 8, this might be due to that some errors in the parameters have 
opposite effects which leads to a kind of improvements in the rate of allocation 
when taking all the errors in the parameters in one case. 
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[4.4] Inversion Testing for Tw -Phase Tw Layers Producing Wells 

The target of this part is to test the inversion for the case of two phase two layers 
production to see whether it is possible to allocate the rate under perfect knowledge of 
input parameters or not The scenarios that will be studied are scenarios 7 and 8 with and 
without imposing the total rate for each phase as constraint in the objective function. The 
true temperature will be generated similarly as have been done before using the true rates 
from the two producing zones to be as follows: q 0 i = 300 STB/day, qo2 = 200 STB/day, 
q w i = 350 STB/day, and q w2 = 300 STB/day. Also the time and the separation between 
the zones will be the same as that presented in scenarios 7 and 8, while the other input 
parameters will be given in Table 4.3. 



Table 4.3: Input data for the synthetic examples of the two phase two layers 
production well for scenarios 7and 8 m 





Value 


Units 


Formation Thermal Properties 






Geothermal gradient, Gj 


0.0274 


Fffl 


Thermal diffusivity, a 


0.04 


ft A 2/hr 


Thermal conductivity, K» 


33.6 


BTU/D-ft-F 


Bottom hole temperature, Tan* 


237.2 


F 


Fluid properties 






Oil API 


30 


API 


Oil specific heat capacity, Cpo 


0.485 


BTU/lbm-F 


Oil viscosity, po 


1.06744 


cp 


Water specific gravity, f* 


1 


BTU/lbm-F 


Water specific heat capacity, Cp* 


1 


Water viscosity 


0.31 


cp 


Well Data 






Total length of the well, L (Depth of the lower zone) 


6792 


ft 


Inside radius of the tubing, r fl 


0.9075 


in 


Outside radius of tubing, r to 


1.1875 


in 


inside radius of the casing, r d 


2.506 


in 


Outside radius of the casing, r M 


2.75 


in 


Wen bore radius, u 


3.75 


in 


Annulus water thermal conductivity, K^*, 


9.192 


BTU/D-ft-F 


Cementing material thermal conductivity, kcem 


96.5 


BTU/D-ft-F 


Tubing roughness 


0.001 


ft | 


Inclination angle to the horizontal, 8 


90 


Deg 
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Table 4.4 and Fig. 43 show the result of the inversion for the two scenarios starting by 

very small values of the rates as an initial guess and taking into consideration that we are 

i 

using only the points that are close to the producing interval in the optimization. 



Table 4 J: Test of the inversion for two-phase two layers production for high and low 
temperature contrasts between the zones with and without imposing total phase rates as 
constraint 







q<>2 




4w2 


qoi 


4o2 


qwi 


Qw2 




Error, 


Error, 


Error, 


Error, 


Error, 


Error, 


Error, 


Error, 




% 


% 


% 


% 


% 


% 


% 


% 




High Temperature Contrast 


Low Temperature Contrast 


Without 


45.2 


24.3 


16.9 


6 


44.75 


30.31 


15.65 


6.6 


Total rate 


















With. 


0.4 


0.6 


0.17 


0.2 


10.2 


15.3 


4.25 


4.96 


total rate 
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Perfect Knowledge of Input parameters 




Without total Rate as 
constraint 



With total rate as 
constraint 



m Error Qo1J*figh 
contrast 

■ Error Qo2 J-figh 
contrast 

a error Qw1_high 
contrast 

■ Error Qw2_high 
contrast 

■ error Qo1 Jow 
contrast 

■ error Qo2Jow 
contrast 

■ error Qwi Jow 
contrast 

■ error Qw2 Jow 
contrast 



Fig. 4.3: Error in rate allocation for high and low temperature contrasts between the 
zones under with and without imposing total rate as constraint for a perfect knowledge of 
input parameters 
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From the results of the above study, it was shown that imposing the total rate of each 
produced phase improves the accuracy of the rate allocation for the high temperature 
contrast but it does not do a good job for the low one, however we assumed a perfect 
knowledge of the input data. That is why it is very difficult to allocate the rates for the 
two phase production for imperfect knowledge of the input data as to accurately allocate 
the rates, any extra information other than total rate for each phase is required to make 
the problem well constraint and not to be trapped in any local minima as was done in the 
above synthetic examples for the low temperature contrast with and without total phase 
rate as constraint and in the high temperature contrast without total phase rate constraint. 
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Chapter 5 
FIELD EXAMPLE 



The primary objective of this chapter is to test the forward modeling with respect to real 
field data and the secondary objective is to see how the inversion will work in handling 
true data other than the synthetic examples that was mentioned in the previous chapter 
and also to check whether the findings from the synthetic examples are reasonable or not 
So, this chapter is organized as follows: the first part will concern with the data 
description and the second part is data analysis, which is devoted to analyze the data and 
discuss the results. 
[5.1] Data Description 

The field example is taken from Bugpn-07 well, which is a 4 zone intelligent well 
completion with 3 ICV's to allow any or all of the producing intervals to be flowed at the 
same time, the completion for this well is given in Fig. 5.1. In March 13, 2003 well 
testing has been done to test the flow from the middle two zones while the other two 
zones are closed, the test shows a total oil production of 769 m 3 /day (4836.5056 bbl/day), 
and total water production of 236 m 3 /day (1484.2852 bbl/day). Fiber optic line has been 
installed in the completion to measure the temperature inside the well bore and also 
information about the geothennal gradient has been given for this area. Fig. 5.2 shows the 
temperature recording inside the well bore, the geothennal temperature for the section 
that is closed to the two producing zones as was mentioned before that we are only 
interested on the data very close to the producing intervals and also the perforation 
interval for each zones. As seen from Fig. 5.2 that the temperature data are very rough, so 
I smoothed the data using a moving average statistical technique every 10 measurement 
points and the smoothed temperature is also plotted in Fig. 5.2, the smoothed data will 
then be used for our analysis. It should be noted that Bugan 07 well is a deviated well and 
the average inclination angle with horizontal for the section that we are interested on is 
46.36°. From the geothennal gradient, the inclination angle, and the bottom hole 
measured temperature, it was found that the temperature contrast between the two zones 
is about 7 F, which is an intermediate case compared to the synthetic examples that have 
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been showed before, where we showed cases of 13 F and 3 F temperature contrast 
between the producing zones. 



8330 



9B2S 



6.180 riwhang* 



m ISV Casing Shoe 



2S1CBQB 




Fig. 5.1: Completion ofBugan 07 Well 
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Temp v F 




Fig. 5.2: Temperature measurement inside the well bore, geothermal temperature, 
smoothed temperature, and the perforation intervals for each producing zones 
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[5.2] Data Analysis 

From a quick look qualitative interpretation of the temperature data whether the 
smoothed or the rough one, it could be easily noticed that not all the perforated interval 
from the upper zone is produced as we can see that only small portion of the upper 
producing interval has a decline in temperature due to the mixing of the fluid of low 
temperature from the upper zone and high temperature from the lower zone which cause 
mis reduction in temperature, while the other portion of the zone does not show mis 
reduction. 

Table 5.1 shows the data used for the forward modeling of temperature of mis field case. 
The data in bold are the data that there are no information about them, so I add those data 
as parameters in the inversion to be changed to match the measured temperature. Values 
of the data given in the table are the values after inversion. For the inverse modeling to 
get the rates from each zone, I used the same objective functions given in Eq. 4.1, where 
the vector (m) are the rates from each zones and the input data given in Table 5.1, and the 
constraints are given in Table 5.2. 

It should be mentioned that the main difference between the input data for the synthetic 
cases, Table 4.3, and for the field case, Table 5.1, is that in field case, there are number of 
temperature measurements inside each producing zone, so the information about those 
numbers should be given, while in case of synthetic examples, we assume that there is no 
measurement points inside the producing zones, so no additional information about the 
number of temperature measurements are required. Jh addition, Table 5.1 needs 
information like the formation density and specific heat capacity to calculate the 
formation diffusivity, while in the synthetic examples; we just input the formation 
diffusivity, just to simplify the calculation. 

Also, the bounded constraints given in Table 52 are taken from the literature 4, 15 so that 
there is no unrealistic results will be obtained to the parameters entered in the 
optimization as a result of inversion. 
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Table 5 J: Input Data for the field example 



brmation Thermal Properties 



Geothermal gradient, 6t 
Formation Density, p« 
Formation heat capacity, C« 
Thermal conductivity, K, 

Bottom hole temperature, T^t h 

Fluid properties 



Oil API 

Oil specific heat capacity, Cpc 

Oil viscosity, u« 

Water specific gravity, 

Water specific heat capacity, Cp* 

Water viscosity 

Well Data 

Total length of the well. L (Depth of the lower zone) 

Inclination to the horizontal, 6 

Inside radius of the tubing, r 8 

Outside radius of tubing, r to 

Inside radius of the casing, r d 

Outside radius of the casing, fa 

Well bore radius, r^ 

Annulus water thermal conductivity, Kan,, 
Cementing material thermal conductivity, kcm 
Tubing roughness 

No. of temperature measurements in the tower producing zone 
No. of temperature measurements in the upper producing zone 
Depth to the top of the upper producing zone 
Depth to the bottom of the upper producing zone 
Depth to top of the top of the lower producing zone 
Production time 



Value 



0.01538 
176.85 

0.7 

96 
191.91 



30 
0.4 

1.06744 
1 

0.477 

0.31 



14089.46 
46.36 

1.5 

1.75 

2.75 

3.5 

4.5 

10 
9423 
0.001 

90 

12 
13428 
13464 
13795 

840 



Units 



F/ft 
Lbm/ft 3 
BTU/ibm-F 
BTU/D-ft-F 
F 



API 
BTU/lbm-F 

CP 

BTU/lbm-F 
cp 



ft 
Deg. 
in 
in 
in 
in 
in 

BTU/D-ft-F 
BTU/D-ft-F 

ft 
Meas. 



ft 
ft 
ft 
days 
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Table 5.2: The constraints for the model parameters 





nan ualifv 


Valus 


Units 


MQl 




0 


STB/D 




>s 


n 

w 


STB/D 




>s 


0 


STB/O 




>s 


A 
\j 


STB/D 


Pe - 


>?5 


75 


Ibm/ft3 


Pe.. 


<S5 


185 


lbm/ft4 


c 


>— 


0.1 


BTU/lbm-F 


p 


<— 


07 


BTU/lbm-F 




>=s 


10 


BTU/D-ft-F 


V. 

IC 


<= 


96 


BTU/D-ft-F 




>= 


8 


BTU/D-ft-F 


Kanw 


<= 


10 


BTU/D-ft-F 


c w 


>= 


0.4 


BTU/lbm-F 


c w 


<= 


1 


BTU/lbm-F 


Cow 


>= 


0.4 


BTU/lbm-F 


Cow 


<= 


1 


BTU/lbm-F 


q ft T 




4836.5056 


STB/D 






1484.2852 


STB/D 



Fig. 5.3 shows ihe comparison between the measured and the modeled temperature after 
the inversion. The rate allocation from each zone were found to be as follows: q 0 | = 2850 
bbVday, q o2 = 1986.5 bbl/day, q w i = 1 174.1 1 bbl/day, q w2 = 310.17 bbl/day 



Tf<F) 




Fig 53: Comparison between modeled and measured temperature 
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Fig. 5.4 shows the absolute errors between the modeled and the measured temperature 
across the section of the well under study, 



From Figs. 5.3 and 5.4, it should be noted that the modeled temperature is in good 
agreement with the measured temperature with a maximum difference less than 1 F. As 
seen fiom Fig. 5.4 that the higher error is in between the two producing zones, which 
might be due to the presence of the IC valve that complicates the physics of the problem 
and we did not take that into account in the modeling. Also for the inverse modeling, 
from tests that have been done before on each of two zones separately, it was noticed that 
the gross production from the lower zone is much higher than for the upper zone for both 
oil and water and that was captured from the inversion results. However, it is worth to 
mention that inversion was found to be highly depend upon the initial guess and that the 
problem is non unique problem and this finding support the results obtained from the 
synthetic examples where we showed for the two phase problem and small temperature 
contrast between the two zones it is very difficult to allocate the rate with high 
confidence due to the non uniqueness of the problem. 



1.2 




Measured Depth, ft 



14200 



Fig. 5.4: Error comparison between the measured and the modeled 
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Chapter 6 

CONCLUSIONS & RECOMMENDATIONS 



[6.1] Conclusions 

From this study, it can be concluded that: 

1- A forward model of temperature has been developed from a basic theory of mass, 
momentum, and energy balance equations to describe the temperature profile in a 
single phase single layer producing well and tested against the numerical and the 
well known models in the literature and the results shows its good accuracy close to 
the producing intervals 

2- The forward model has been extended to two zones two phase (oil-water) producing 
wells to be used as a forward tool in inverting for the zonal rate by knowing the 
temperature from DTS measurements 

3- An inverse modeling using the GRG optimization algorithm coded in Excel is used 
to invert for the zonal rate allocation by minimizing the difference between the 
measured temperature from the DTS measurements and the calculated from the 
forward model 

4- Several synthetic examples have been studied to test the validity and the accuracy of 
zonal rate allocation from DTS measurements in two layers producing wells under 
different conditions and the results reveal the following: 

i- For single-phase liquid production with high temperature contrast 
between the producing zones, the zonal rates can be allocated with good 
accuracy without imposing the total rate as constraint in the optimization 
problem 

ii- For single-phase liquid production with low temperature contrast 
between the producing zones, the zonal rates can be allocated with good 
accuracy if the total rate is added as constraint in the optimization 
problem. 

iii- For two-phase (oil-water) production with high temperature contrast 
between the producing zones, the zonal rates can be allocated with good 
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accuracy only after imposing the total rate for each production phase as 
constraint in the optimization problem, 
iv- For two-phase (oil-water) production with low temperature contrast 
between the producing zones, the zonal rates are difficult to be allocated 
even if the total rate for each phase is added as constraint, the problem 
shows a high non-uniqueness and the optimization mainly depends on the 
starting guess of the rates. So, prior information to select a good starting 
guess for the rates, or another constraint can be added to the problem to 
make it well constraint and this might improve the accuracy of the 
inversion. 

5- A simple Excel spreadsheet has been prepared to analyze a field DTS measurements 
from two zones two layers producing wells. 

6- Both the forward and the inversion has been tested on a field data taken from Bugan- 
07 deviated well and the results show a good agreement of the forward model with 
the measured temperature from DTS within maximum difference of 1 F. Also the 
results of inversion show qualitatively a good agreement with the information about 
the zonal productivity. 

[6.2] Recommendations 

The following points are recommended for future work: 

1- The Excel spreadsheet that was developed to handle the field case for two zone 
two phase producing wells can be extended for more than three zones and the 
extension is trivial. This allows for analyzing the DTS measurements from wells 
producing two phase (oil-water) or single phase (oil) from one, two or three zones 
in short time. 

2- Extending the forward modeling for single-phase g?ts and multi phase flow (gas 
and liquid). The extension of the forward model for multi phase flow (gas and 
liquid) can definitely handle a single-phase gas production, thus the key is to 
extend the forward model to multi phase flow. The major difficulties in modeling 
a multi phase arising in finding a way to get the Joule Thomson coefficient and 
the pressure drop inside the well to substitute them in Eq. 2.26. Reference 7 has 
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mentioned a good way to get an approximate equation for the Joule Thomson 
coefficient for black oil models and also for the pressure drop. The pressure drop 
can be calculated from any multi phase pressure drop correlation while the Joule 
Thomson coefficient can be obtained as a rate weighting average between the 
liquid Joule Thomson coefficient obtained from Eq. 2.30 and that for gas which 
depends upon a first order derivative of the compressibility factor with respect to 
temperature, any correlation for the gas compressibility factor can be used to get 
this first order derivative. 
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